Matrix Product State and mean field solutions for one-dimensional systems can be 

found efficiently 
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We consider the problem of approximating ground states of one-dimensional quantum systems 
within the two most common variational ansatzes, namely the mean field ansatz and Matrix Product 
States. We show that both for mean field and for Matrix Product States of fixed bond dimension, 
the optimal solutions can be found in a way which is provably efficient (i.e., scales polynomially) . 
This implies that the corresponding variational methods can be in principle recast in a way which 
scales provably polynomially. Moreover, our findings imply that ground states of one-dimensional 
commuting Hamiltonians can be found efficiently. 



I. INTRODUCTION 

Characterizing the ground states of quantum spin sys- 
tems is a highly challenging task. Different from the sit- 
uation for classical systems, where for the very least, the 
ground state can be described efficiently, this does gen- 
erally not hold for quantum systems. Thus, their ground 
states are considerably harder to compute: as Kitaev has 
shown, solving this problem is likely to be hard even for 
quantum computers [l|, Hj]. Even more surprising, finding 
the ground states of quantum systems remains equally 
hard when we restrict our interest to one-dimensional 
systems [|| . This is in sharp contrast to the case of clas- 
sical spin systems, which can be efficiently solved in one 
dimension, whereas two-dimensional systems are known 
to be NP-hard. 

Despite these hardness results for ground states of 
quantum systems, physical properties of interest can fre- 
quently be determined efficiently using numerical meth- 
ods. Even a mean field ansatz, which completely ne- 
glects correlations between the particles, may already 
reproduce many physical quantities relatively well [J]. 
In most cases however, correlations arc important and 
other methods must be applied. While in two dimen- 
sions, imposing frustration or fermionic statistics yields 
Hamiltonians which cannot be assessed well by numeri- 
cal methods, one-dimensional systems - despite the hard- 
ness results for ID Hamiltonians - generally turn out to 
be extremely well simulatablc using a method known as 
the the Density Matrix Rcnomalization Group (DMRG) 
algorithm [j| DMRG can be understood as a vari- 
ational method over the class of Matrix Product States 
(MPS) @,i : These states can be considered as a gener- 
alization of the mean field product states, but with a lim- 
ited amount of correlations. It has been proven that this 
allows for the efficient approximation of ground states of 
gapped local Hamilonians Q. 

The motivation for this work originates in the following 
observation: On the one side, DMRG is a highly success- 
ful algorithm which finds the correct minimum in basi- 
cally all practical cases. On the other hand, it has been 
shown that the problem of finding MPS ground states can 
be NP-hard [Tl| , as well as certain minimization problems 



encountered in the DMRG algorithm [lOj . The contrast 
between these hardness results and the apparent success 
of mean field and MPS algorithms raises the question 
whether it is possible to prove the efficiency of variational 
methods over the class of mean-field and Matrix Product 
States in full generality for gapped quantum systems, or 
whether even studying physical ID problems is hard, as 
it is the case in two dimensions. 

In this paper, we settle this question by studying 
whether (and to which extent) for a given Hamiltonian 
it is possible to find the optimal mean field state or MPS 
for a fixed bond dimension. Surprisingly, we find that the 
same technique which is used to show that classical spin 
chains can be solved efficiently also allows for solving the 
problem in the case of mean field theories and MPS. 

In particular, we show that it is always possible to 
find the optimal mean field solution of a qu-<i-it chain of 
length N up to accuracy 1/e in energy in a time which 
scales as (N/e) M . Concerning Matrix Product States, 
we find that approximating the optimal MPS of bond 
dimension D up to precision I/e requires a computation 
time which scales as as {N 2d /e) 6D , where D is the bond 
dimension of the MPS. For a fixed bond dimension, the 
scaling is thus polynomial both in the length of the chain 
and the accuracy. 

Note that the result for MPS (in particular, the ex- 
ponential scaling in D) is essentially optimal, since is 
has been shown that the difficulty of the problem has to 
scale exponentially with D (more precisely, it is NP-hard, 
where D is polynomial in the problem size [HI). Note 
also that while the exponential scaling seems daunting, 
in practice one is typically interested in intensive quan- 
tities for which a moderate D will suffice [l2j]; moreover, 
the polynomial scaling of the algorithm in N typically 
allows for the efficient evaluation of such quantities even 
in the thermodynamic limit (2lj . 

Our findings show that variational calculations over 
both the mean field ansatz as well as Matrix Product 
States can be carried out in a way which is promised 
to yield the optimal solution, thus resolving the ques- 
tion as to whether variational algorithms over MPS - 
such as DMRG - can be rephrased in a way which prov- 
ably succeeds. Moreover, while the algorithm might be 
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unpractical for D's of several 100 as used in practical 
DMRG implementations, the algorithm could prove use- 
ful in practice to find the optimal MPS ansatz for a low 
D, which can then be used to bootstrap DMRG [lH, thus 
helping to avoid local minima. Finally, our findings also 
imply that ground states of one-dimensional commuting 
Hamitonians can be found efficiently, as they can be ex- 
pressed as MPS. 



II. CLASSICAL SPIN CHAINS 

We start by explaining how to solve a classical spin 
chain efficiently While this is known, it will help us to 
clarify the essential ideas of the proof technique, as used 
later for mean field and MPS. Given an open bound- 
ary condition (OBC) Hamiltonian hk t k+i(ik, *fe+i), ik £ 
{1, . . . , d}, k = 1, . . . ,N, we want to find the i\, . . . ,ijv 
which minimize the energy, 

E = min [h 12 (ii, i 2 ) H h /ijv-i,Jv(*JV-1,«Jv)] ■ (1) 

tl,...,«JV 

To this end, define recursively 

E 1 {i 1 ) = (2) 

E k (ik) = min [B fc _i(i fe _i) + h k -i,k{ik-i, ik)} , k>2. 
ik-i 

Then, the ground state energy is given by E = 
minj N En(in), and the minimization can be carried out 
efficiently by computing the Ek{ik) of Eq. ((2]) sequen- 
tially: The reason is that in every step, Ek (the minimal 
energy of the half chain left of k) only depends on the 
value ik of the spin at k - this is the only variable which 
we still need to access to minimize the energy of Hamilto- 
nian terms to come. The computational cost is as follows: 
for each of the d values of ik, one has to compute Ek{ik)- 
Each computation involves the minimization over d set- 
tings of ifc-i, and the total computational cost is thus 
Nd 2 . Note that the algorithm not only yields the op- 
timal energy, but also the corresponding ground state 
ii, ... ,in- 

The intuition behind the algorithm is to proceed from 
left to right through the chain and at every site minimize 
the energy of the half chain left of k as a function of the 
boundary setting. Here, the "boundary setting" contains 
all those spins whose value will still influence the optimal 
energy of Hamiltonian terms to the right of k; in the 
case of the classical system, this is just the spin at the 
boundary. Then, the optimization can be carried out 
sequentially by adding one new interaction at a step and 
minimizing the total energy of the left half chain (i.e., 
the previous total energy plus the new coupling) as a 
function of the new boundary. Sloppily speaking, the 
idea is that while proceeding through the chain, we have 
to make choices about the state, and we want to keep 
the dependence of the minimal energy of the half-chain 
on all past choices which can still influence the future. 



This construction contains all important ideas for the 
mean field as well as for the MPS setting. For mean 
field, the boundary condition is again only the spin on 
the boundary, which however is now a continuous degree 
of freedom. Thus, we have to show that we can discretize 
its values without loosing too much accuracy. For the 
case of MPS, the situation seems more involved, since all 
choices in the past can influence the future. However, 
as we will show later, MPS have exactly the property 
that the influence of past choices on the future is fully 
characterized by what is passed through the bonds, and 
is thus bounded. 



III. MEAN FIELD 

Having understood the method, let us now consider 
the problem of minimizing the energy of a ID quantum 

system H = J2k=i N-i H k,k+l, \\H k ,k+i\\oo < 1 with 

respect to a mean field ansatz \ip) = ® fe=1 \ipk), IV'fc) £ 
C d (in the following, all states are normalized) , where we 
try to minimize 

2V-1 

E= min (ip k , ipk+i \ H k ,k+i \ipk, i>k+\) ■ (3) 

This minimization again allows for a recursive formula- 
tion as in @. However, since the parameters \ipk) arc 
continuous, and the cost functions Ek(\ipk)) are non- 
linear and thus cannot be solved for exactly, we restrict 
the minimization to an e-net, i.e., a discrete set of \ipk)i 
a = 1 , . . . , A such that 

V|Vit) eC 3a: ||V*-V*l|i<e 

(we use the convention <fi = \<fi) (4>\). As shown in [l4j 
(Lemma II. 4), such a net of size A < (5/e) 2d exists. This 
reduces the algorithm to the algorithm for classical ID 
chains described above, which will yield the optimal so- 
lution in the set of all product states (^) IV'fe) on the net. 

The proper way to think of the net is as a constraint 
made on the mean field ansatz as a whole, rather than 
of a way to approximate each recursion step separately 
(since this would lead to accumulating errors). Thus, in 
order to bound the error made by introducing the net, we 
just have to bound the difference between the minimum 
in the set of product states and the minimum on the net 
in Eq. ([3]). To this end, start from the optimal product 
state §Z)\ipk) an d replace each state by an e-close state 
\ipk) on the net. For each term in the Hamiltonian, this 
yields an error of at most 

\tr[H k , k+l (^ k 9 - r k ik) ® ^| fc ! +1) )]| 

< \\H kt k+i\\oo\\i>k®i>k+i-^ {k) 

< 2e . 

Thus, the total error from restricting the minimization of 
the energy © to states in the net is 2Ne. For a targeted 



3 



accuracy 8 in energy, we thus have to choose e = 8/2N, 
so that the size of each net will be A = (lON/8) 2d . The 
minimization can be thus rewritten in the iterative form 
© and carried out in time d 4 A 2 N = Nd 4 (10N/8) 4d . 



IV. MATRIX PRODUCT STATES 

Let us now turn towards matrix product states. 
Starting again from a ID quantum Hamiltonian H = 

J2k=i Hk,k+i, \\Hk,k+i\\oo < 1, we wish to minimize the 
energy over all Matrix Product States [l5[ 



\x({A k })) = tr[Al...A»]\i 



of a given bond dimension D, A k G Mdxd (except for 
A\ G M\xD, Af G Mdxi)- It is known (cf. that 
every MPS can be brought to a standard form for which 
J2i A i( A iY = 1- With this gauge, |x({-4 fc })> is normal- 
ized, and the energy of a term Hk t k+i can be computed 
as follows: Define 



Pi = l, 

Pfc+1 =7e(A fc ,p,):-^(Af)W« fe ■ 

Then, the energy of Hk 3 k+i is given by 

£k(A k ,A k+1 , P k) = 

J2 (a,b\H ktk+1 \c,d)tT[(A k a + y(A k b y Pk A k c A^} 

a,b,c,d 



(4) 



(5) 



The task is to minimize the total energy over the set of 
MPS, 



E 



N-l 

ruin > 

A 1 ,....!" ^ 

' ' fe=l 



£ fc (A*,,4 fc+ \p fc ) 



(6) 



Note that by virtue of the recursion relation Q, ac- 
tually depends on all A 1 with I < k + 1, thus seemingly 
ruling out the same approach as for the mean field ansatz. 

To resolve this, we rewrite the minimization over all 
Ak in ^ as a series of constrained minimizations, 



and pk at the boundary, and thus solve the minimization 
problem ©. 

Clearly, it will again be necessary to use nets to be able 
to implement the optimization efficiently. We will put a 
net on both the A's and the p's, and denote elements of 
the nets by A and p. Let us define an operation Af : 
p i y p which maps every p to the closest point in the 
net. We define a notified recursion relation for the p's, 
1Z = AfolZ, and define all minimizations (in particular the 
constrained minimizations Ek) with respect to variables 
A k and pk in the net, and constraints according to the 
notified recursion relation 1Z • 

This coarse-grained version of the iterated protocol, 
which can be carried out efficiently, will give the opti- 
mal solution in an ansatz class which is defined by the 
coarse-grained variables and the coarse-grained recursion 
relation - note that this is not equal to the optimal solu- 
tion with respect to the coarse-grained MPS due to the 
additional coarse-graining Af in the the p's. 

To bound the error made by introducing the nets, 
we compute how much the energy of an arbitrary MPS 
changes due to the coarse-graining of the ^4's and the 
p's. To this end, we first bound the difference in energy 
caused by coarse-graining the p's as compared to the en- 
ergy of the MPS described by the same A k 's, and second, 
we bound the error in energy made by coarse-graining the 
A k, s. By choosing the nets such that both of these errors 
become small, we make sure that the A k 's found by opti- 
mizing the above recursion relations yield an MPS which 
is close in energy to the optimal MPS. 

Let us first bound the error made by inserting Af. We 
will put an e p -net on the p's, i.e. for each p there is a p 
in the net with ||p — p||i < e p . From (H}, 

\\pk - Pk\\i = 

= \\Af(j2( A i) f h-iAf) -^(iftWi^lL 
= l^^iPk-i-Pk-i^ + ^ll (Hi < i) 

<\\j2A k (A k )Hp k -i-Pk-i: 

< \\Pk-i - Pfc-i Hi + e p < • • • 

< (k - l)e p . 



mm = mm 

Ai,...,A N An-iAn,Pn- 



mm 

(Ak-l,Pk-l)->Pk 



mm 

(A 1 ,p 1 =l)^p 2 



where (Afc_i,pfc_i) -> pk denotes the tuples (Ak-i, Pfe-i) 
which arc compatible with pk, i.e., 1Z(Ak—i, Pfe-i) = Pk- 
Based on this rephrased form of the minimization ([5]) , we 
define 



Ek(A k+ \p k+1 ) 



min £ k (A k ,A' : 

(Ak,Pk)^>-Pk+i 



fe+1 „k 



p k ) + E k -i{A\p k ) 



(if the minimum is over an empty set, let E^ = +oo), 
with Eq = 0. Then, we can sequentially solve for the E k , 
always only keeping track of them as a function of A k 



To bound the effect of the error in p k on each Hamiltonian 
term in Eq. ([5]), note that ([S]) can be rewritten as 

£k(A k ,A k+ \p k ) = tr[V(Hk,k + i<E>^pk} 

with an isometry V a , a bp = (A k A k+1 ) al3 , 

H a ,b,l3 Va,abpV a i^ a bi3 = 8 a ,a' ■ Then, 

\\£k(A k ,A k+1 ,pk)-£k( Ak , Ak+ \pk)\\i = 

= \\tr[V(H k ,k+i®WHpk-pk)]\\i 

< \\V(H k , k+1 ® l)yt|| 

oo||Pfc — Pfc||l 

< (k - l)e p . 
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Thus, the total error in energy due to the net on the p's 
is 



N-i 

E 

k=l 



(k - l)e p < \N 2 e p 



Second, we have to bound the error in energy made 
by replacing the A k by A k chosen from an e^-net, which 
approximates the A k in operator norm up to €a- To this 
end, define A fc = A k - A k and E k (p) = Y, l A k p{A k )\ 
= E t A k p(A k y. Then, E k (l) = 1, E k is contrac- 
tive with respect to || • Hoc, and 



|Dfe(p)llc 



< 



Mp\U\&i\\oc<de A \\p\\ 



The overlap of the MPS with and without a net is given 
by 

(x({A k })\ X {A k }) = ((Ex + ©i) o • ■ • o (E N + D A r))(l) . 

For the iterated application of (E k +Bfe), we use that the 
deviation a k from the identity grows according to 



kfc-i||oo = \\(E k +B k )(l+a k ) - t\ 

< ||E fc (tr fc )||oo + ||Bfc(l)||oo 

< (1 + d£A)\\<rk\\oo + &A ■ 



Under the condition that ||tTAr— fc ||co < 2kdtA < 1, this 
yields \(x{{A k })\ X {A k }}\ > 1 - 2Nde A - From this, we 
can derive a bound on the difference in energy, 

6A = \tr{H( X ({A k })- X ({A k }))}\ 
<\H\\ 00 \\ X ({A''})-xaA''})\\ l 



< 2N\ 1 - 



< 4iV 



(x({A k })\ X {A k }) 



Let the targeted error in energy now be 5, and choose 



So 



8/2. Thus, we will need nets of precision 



8/N 2 and ca = 5 /6AN 3 d, respectively. Such nets 



6a 
e P = 

exist of size A p = (3/e p ) D ' and A A = (3/e A ) 2d£H ([3, 
see Appendix). Each evaluation of an energy in ([5]) takes 
d D 3 elementary steps, and thus the total computational 
cost is 



Nd 4 D J 



g2d+l2l2d^y6d+2^2d 



s 3 



2D Z 



(7) 



V. COMMUTING HAMILTONIANS 

Our results also imply that the ground state of a lo- 
cal Hamiltonian on a one-dimensional chain with mu- 
tually commuting terms can be found efficiently. This 
follows straight away from the fact that the ground state 
of any commuting Hamiltonian can be expressed as an 
MPS with bond dimension D = d 2 fnl. 



Alternatively, one can map the problem of solving 
any one-dimensional commuting Hamiltonian to solving 
a classical ID chain, as shown by Bravyi and Vyalyi [l8[ , 
which allows to use the classical ID algorithm described 
at the beginning. (In fact, they show that every com- 
muting 2-local Hamiltonian - which is always the case in 
ID - can be mapped to a classical problem on the same 
interaction graph.) 



VI. PERIODIC BOUNDARY CONDITIONS 

Up to now, we have focused on Hamiltonians on OBC 
ID chains. Let us now briefly discuss how to adapt our 
method to the case of a periodic boundary condition 
(PBC) ID system. 

In the case of both classical chains and the mean field 
ansatz, the PBC case can be tackled by additionally let- 
ting all E k in ([2]) depend on i\ (\if>?)): the reason is that 
i\ is also part of the half-chain boundary, as the value of 
i\ will influence the optimal energy when evaluating the 
Hamiltonian term ft.jv,i- 

To solve a PBC Hamiltonian with (OBC) MPS, the sit- 
uation is a bit more subtle, as the energy of the coupling 
Hn,i does not only depend on A , but also on the way in 
which it is passed through the chain. Two possibilities to 
deal with that would either be to evaluate the E k not only 
as a function of A k and p k , but also of the p^f arising 
from putting an operator basis \a) (/3| at site 1 [i.e., start- 
ing the recursion from the state p% = {A 1 a )' 1 Ap]; or to 

keep the dependence on A 1 and A k , and instead of keep- 
ing the dependence on p k rather make E k depend on the 
possible channels TZ(A k , TZ(A k ~ 1 , ■■■ , K(A 2 , •))) passing 
p\ through the MPS. A more efficient way to deal with 
PBC Hamiltonians is however to fold the chain: This 
yields an OBC Hamiltonian on a chain of length N ~/2 
which can be solved with the original algorithm using an 
OBC MPS with bond dimension D 2 , which includes the 
case of the folded MPS of dimension D. Similarly, the 
case of PBC MPS with bond dimension D can be solved 
by embedding it in an OBC MPS with bond dimension 
D 2 . 
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Nets 

We use [lH, Lemma 4.10. It states that in K n with any 
norm ||| ■ |||, there exists an e-net with respect to ||| • | 
for the unit sphere with cardinality at most (1 + 2/e)" < 
(3/e)». 

This Theorem can be applied directly to the case of 
D x D density operators which span a Z? 2 -dimcnsional 
real vector space. The size of an e p -net with respect to 
the trace norm is thus A p < (3/e p ) D . To create such a 
net, one can e.g. place a lattice on M. D and use that the 



trace norm is bounded relative to the Euclidean distance. 

For the case of the space of isometrics A : C D — > 
C d x C D , it is necessary to embed them in a vector space, 
namely the space of all complex D x dD matrices. Then, 
the theorem states that there is a net of isometries of 
size at most Aa < (3/eA) 2dD ■ (Assuming that the net 
consists of unitaries imposes no restriction, since the the- 
orem bounds the size of any maximal set of vectors with 
distance > e.) By utilizing the Solovay-Kitaev theorem, 
we can generate such a net efficiently for any D [l], [lj| . 
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